Investigation of the mechanism of Prunella vulgaris in treatment of papillary thyroid carcinoma based on network pharmacology integrated molecular docking and experimental verification

To analyze the molecular mechanism of Prunella vulgaris L. (PV) in the treatment of papillary thyroid carcinoma (PTC) by using network pharmacology combined with molecular docking verification. Traditional Chinese Medicine Systems Pharmacology Database and Analysis Platform database was used to predict the main active components of PV, Traditional Chinese Medicine Systems Pharmacology Database and Analysis Platform, PubChem, and Swiss Target Prediction databases were used to obtain the corresponding targets of all active components. Targets collected for PTC treatment through Gene Cards, Digest and Online Mendelian Inheritance in Man databases respectively. The Search Tool for the Retrieval of Interaction Gene/Protein database was used to obtain the interaction information between proteins, and the topology analysis and visualization were carried out through Cytoscape 3.7.2 software (https://cytoscape.org/). The R package cluster profiler was used for gene ontology and Kyoto encyclopedia of genes and genomes analysis. The “active ingredient-target-disease” network was constructed by using Cyto scape 3.7.2, and topological analysis was carried out to obtain the core compound. The molecular docking was processed by using Discovery Studio 2019 software, and the core target and active ingredient were verified. The inhibition rate was detected by CCK8 method. Western blot was used to detect the expression levels of kaempferol anti-PTC related pathway proteins. A total of 11 components and 83 corresponding targets in the component target network of PV, of which 6 were the core targets of PV in the treatment of PTC. It was showed that quercetin, luteolin, beta (β)-sitosterol, kaempferol may be the core components of PV in the treatment of PTC. vascular endothelial growth factor A, tumor protein p53, transcription factor AP-1, prostaglandin endoperoxidase 2, interleukin 6, and IL-1B may be important targets for the treatment of PTC. The main biological processes mainly including response to nutrient levels, response to xenobiotic stimulus, response to extracellular stimulus, external side of plasma membrane, membrane raft, membrane microdomain, serine hydrolase activity, serine-type endopeptidase activity, antioxidant activity, etc IL-17 signaling pathway, and PI3K-Akt signaling pathway may affect the recurrence and metastasis of PTC. Kaempferol may significantly reduce the activity of Papillary cells of human thyroid carcinoma bcpap cell lines cells compared with quercetin, luteolin, β-sitosterol. Kaempferol may reduce the protein expression levels of interleukin 6, vascular endothelial growth factor A, transcription factor AP-1, tumor protein p53, 1L-1B and prostaglandin endoperoxidase 2, respectively. PV has the characteristics of multi-components, multi-targets and multi- pathways in the treatment of PTC, which network pharmacology help to provides a theoretical basis for the screening of effective components of PV and further research.


Introduction
Globally, among the malignant tumors of the thyroid endocrine system, thyroid cancer has a high incidence rate, accounting for more than 80% to 90% of all endocrine tumors. [1] Among thyroid cancers, papillary thyroid cancer (PTC) has the highest incidence rate, accounting for more than 80% of thyroid malignant tumors. The formation of PTC is a very complex biochemical process, characterized by the abnormality of various molecules. [2] Large scale genome research shows that genetic changes play a key role in the occurrence and development of PTC, but the reasons for its continued high incidence are not yet clear. [3] Thyroid disease screening has become a routine examination item of health examination, but there is still a lack of thyroid disease screening in some people's health examination items, and some people do not have the awareness of regular health examination, and they do not go to the hospital until they find a neck lump or a symptom of neck discomfort. [4] At present, total thyroidectomy, radioiodine therapy and thyrotropin suppression therapy have been highly developed and widely accepted. [5] The adverse reaction of Traditional Chinese Medicine (TCM) is relatively small, which can help reducing the recurrence and metastasis of tumors, improving the quality of life of patients, with obvious advantages. [6] Prunella vulgaris L. (PV), a traditional Chinese herbal medicine in China, is a plant of Labiatae. It is distributed around the world, mainly in temperate and tropical mountains, such as Europe, Asia, northwest Africa and North America. It is mainly distributed in the Huaihe River basin and the middle and lower reaches of the Yangtze River basin in China. [7] According to the Pharmacopoeia of the People's Republic of China, it has the effects of clearing the liver, purging fire, improving eyesight, dispersing knots, and reducing swelling in clinical efficacy. [8] Therefore, PV is also commonly used in clinical treatment of thyroid cancer, breast cancer, liver injury, hypertension, hyperglycemia, prostate diseases, and other diseases. [9][10][11][12] PV mainly contains terpenoids, phenolic acids, organic acids, sterols, coumarins, flavonoids, polysaccharides, and volatile oils. [10] Modern pharmacological studies have shown that its extracts have a variety of biological activities, such as antioxidation, regulating immune function, anti-inflammatory, antibacterial, antiviral, and regulating the microenvironment of tumor metastasis. [10,11] With the extensive clinical application of PV, the experimental research of PV has gradually deepened, and its pharmacological effects and molecular mechanism of action are gradually clear. However, the active ingredients contained in PV are relatively complex, and its main medicinal substances and pharmacological mechanism are still unclear. The inhibition of thyroid cancer, active ingredients and potential action targets of PV are still lack of systematic research.
Network pharmacology can successfully establish the relationship among "drugs-components-targets-diseases" and reveal the complex functions of Chinese medicine at multiple levels and systems. [13] The network pharmacology technology can be used to obtain the relationship between "multi components and multi targets" and protein interaction, predict the TCM and its molecular mechanism of action, and make the research of TCM extend better. This study adopts network pharmacology combined molecular docking methods to explore the effective components and potential mechanism of action of PV in inhibiting PTC, so as to provide more theoretical support for clinical use of PV in treating PTC.

Prediction of main active components and corresponding targets of PV
Using Traditional Chinese Medicine Systems Pharmacology Database and Analysis Platform (TCMSP) (http://tcmspw. com/tcmsp.php) to search the chemical components of PV, select the components with oral bioavailability≥30% and drug like property≥0.18 as the effective components, download the corresponding chemical structure, and save it in MOL2 format.

Construction of the relationship of active ingredient and disease target
Through TCMSP database, Swiss Target Prediction (http:// www.swisstargetprediction.ch/) database and PubChem database (https://pubchem.ncbi.nlm.nih.gov/) respectively. The target library corresponding to the main active components of PV was obtained by searching the corresponding targets of all active components in PV and removing duplicate genes. Using the Gene Cards database (https://www.genecards.org/), Disgenet database (https://www.disgenet.org/), and Online Mendelian Inheritance in Man database (https://omim.org/) respectively to collects targets for PTC treatment with the key word "papillory thyroid cancer" OR "PTC," removes duplicate genes, and obtains the corresponding target library.

Protein-protein interaction (PPI) network construction
In order to further study the action mechanism of PV action target and disease target at the protein level, we use R language to obtain the intersection between them, and use R package Venn Diagram for visualization. Input the intersection gene into Search Tool for the Retrieval of Interaction Gene/Protein database (https://string-db.org/) to obtain the interaction information between proteins, and conduct topological analysis and visualization through R language.

Gene ontology (GO) and Kyoto encyclopedia of genes and genomes (KEGG) analysis
In order to further study the core target function and the main action pathway of PV in the treatment of thyroid cancer, the R package cluster profiler was used to conduct GO and KEGG analysis on the intersection targets. The threshold values were set at P<.05, q<0.05, and finally 20 results with the lowest P value were retained for visualization. To validated the anti-PTC mechanism of PV across the key targets and multiple pathways, the KEGG mapper functional analysis was utilized to mark the target genes on the pathway associated with PTC.

Build a network of "active ingredient -target-disease"
In order to analyze the correlation between active compounds, a network of "active ingredient-target-disease" was constructed by using Cytoscape 3.7.2 software (https://cytoscape.org/). The function of "network analyzer" is used for topological analysis to obtain the core active ingredients.

Molecular docking
In order to ensure the reliability of molecular docking, the protein crystal structure with resolution < 2.5 Å and ligand complex is selected to establish a molecular docking model. The molecular docking program uses the Libdock module of the Discovery Studio 2019 software package. As a small molecule ligand, the key active ingredient CAS number is searched through the PubChem database, derived from Chem 3-dimensional 19.0, and the lowest energy is calculated and saved in MOL2 format. As the receptor, the core target protein is downloaded from the protein data bank (PDB) database in PDB format. First, the water molecules in the PDB structure are removed, and the docking active pocket is defined by the original ligand molecules. After setting the docking parameters, the ligand molecules in the crystal structure are extracted and re docked to the predefined active pocket. At the same time, the root mean square difference between the docked ligand molecule conformation and the initial conformation in the crystal structure is calculated. If the root means square difference value is <2.5, the molecular docking result is considered reliable. On this basis, molecular docking was performed on the core target and corresponding compounds to verify the binding strength of the core target and the main active components of PV.

Statistical methods
The data was analyzed by SPSS 20.0 statistical software (SPSS Inc., Chicago, IL), and all the images were drawn by GraphPad Prism 9.0 software (https://www.graphpad.com/). Image J software (https:// imagej.net/software/imagej/) was used for protein gray scale calculation, and P<.05 was considered as statistically significant difference.

Main active ingredients of PV
According to the screening conditions (oral bioavailabil-ity≥30%, drug like property≥0.18), 11 main active ingredients of PV were screened (Table 1).

Construction of active ingredient and disease target library
121 active ingredient targets of PV were retrieved through TCMSP, Swiss Target Prediction and PubChem database, and these targets were converted into gene names through UniProt database. Through searching Gene Cards and Online Mendelian Inheritance in Man databases, 3064 PTC related targets were obtained. The R package Venn Diagram was used to finally visualize 83 intersections of active ingredient targets and disease targets (Fig. 1A).

Network diagram of "active components-drugs-target-diseases"
According to the above results, the network of "active components-drugs-target-diseases" was constructed by using Cytoscape 3.7.2 (Fig. 1B). The triangle represents the PV, the circle represents the target, and the hexagon represents the effective component of the PV monomer compound. Nerwork analyzer analysis shows that quercetin, luteolin, beta (β)-sitosterol and kaempferol have a degree value of the first 4 drugs, which is significantly higher than other effective ingredients (Table 2). At the same time, the network diagram shows that multiple targets are related to different active ingredients, which reflects the characteristics of TCM in treating "multiple-components and multiple-targets."

PPI network construction
Input 83 intersection targets into Search Tool for the Retrieval of Interaction Gene/Protein database to get the interaction information among targets, and use the software Cytoscape 3.7.2 for topological analysis and visualization to get the interaction network diagram of genes (Fig. 1C). The ordinate represents the target gene name, and the abscissa represents the connectivity. The core network was screened according to the first 6 principles of degree, and the interaction network diagram of 6 core genes was obtained (Fig. 1D).

Result of GO and KEGG analysis
GO ( Fig. 2A) and KEGG (Fig. 2B) analysis of intersection targets are conducted through R package cluster profiler. A total of 1665 GO entries (q<0.05) were obtained from GO function enrichment analysis, including 1503 biological process entries, mainly including response to nutrient levels, response to xenobiotic stimulus, response to extracellular stimulus, response to oxidative stress, cellular response to chemical stress, response to reactive oxygen species, cellular response to oxidative stress, response to lipopolysaccharide, reactive oxygen species metabolic process, cellular response to reactive oxygen species etc. A total of 49 cellular composition entries mainly including external side of plasma membrane, membrane raft, membrane microdomain, vesicle lumen, ficolin-1-rich granule lumen, ficolin-1-rich granule, secretory granule lumen, cytoplasmic vesicle lumen, plasma membrane raft, caveola etc. One hundred-thirty molecular function entries mainly including serine hydrolase activity, serine-type endopeptidase activity, antioxidant activity, nuclear receptor activity, ligand-activated transcription factor activity, serine-type peptidase activity, RNA polymerase II CTD heptapeptide repeat kinase activity, protein serine/threonine/ tyrosine kinase activity, virus receptor activity, exogenous protein binding etc. By KEGG pathway analysis (q<0.05), 149 pathways were obtained. The first 20 pathways are core pathways, mainly including TFluid shear stress and atherosclerosis, AGE-RAGE signaling pathway in diabetic complications, Lipid and atherosclerosis, IL-17 signaling pathway, Malaria, Prostate cancer, Chemical carcinogenesis-receptor activation, PI3K-Akt signaling pathway, Rheumatoid arthritis, EGFR tyrosine kinase inhibitor resistance, Bladder cancer, TNF signaling pathway, Endocrine resistance, Kaposi sarcoma-associated herpesvirus infection, Breast cancer, Chemical carcinogenesis -reactive oxygen species, Relaxin signaling pathway, Melanoma, Non-small cell lung cancer, Calcium signaling pathway etc. which are mostly related to tumors, indicating that PV may treat PTC by acting on different signal pathways. Annotated map of the key target genes locations of PV in PTC related pathways was presented in Figure 3. It was revealed that most of the key target genes are associate with IL-17 signaling pathway (Fig. 3A) and PI3K-Akt signaling pathway (Fig. 3B) which may affect the recurrence and metastasis of PTC.

Molecular docking results
The 2-dimensional and 3-dimensional structural diagrams of the 4 active compounds of PV were presented in Figure 4A and Figure 4B respectively. Based on the PPI network, we selected the top 6 of degree as the core targets for molecular docking with 4 active compounds of PV (Fig. 4C). It is mainly combined with the target through traditional hydrogen bond, hydrocarbon bond, Pi-Pi conjugation, etc. In addition, the libdock score of each targets and compounds are presented in Figure 4D. Libdock score indicates the binding degree of ligand and target protein crystal. The higher the value is, the stronger the molecular binding ability is. Among them, those with the strongest combination ability with JUN are β-sitosterol (120.529) (Fig. 5A) and kaempferol (104.874) (Fig. 5H), the strongest binding ability with VEGFA are luteolin (71.465) (Fig. 5F), β-sitosterol (52.1039) (Fig. 5C) and kaempferol (71.3777) (Fig. 5I), the strongest binding ability with recombinant human interleukin 1 beta (IL1B) are luteolin (100.753) (Fig. 5D) and kaempferol (80.8784) (Fig. 5J), the strongest binding ability with IL6 is kaempferol (79.3898) (Fig. 5G), the binding ability with PTGS2 are luteolin (88.4373) (Fig. 5E) and quercetin (106.552) (Fig. 5K), the binding ability with TP53 is quercetin (71.63) (Fig. 5L) and β-sitosterol (106.552) (Fig. 5B) respectively. It was revealed that the 4 core components of PV exhibit good binding with PTC related targets respectively (Table 3, Fig. 5).

The IC 50 values of active compounds in PV.
The active compounds quercetin, luteolin, β-sitosterol, and kaempferol from PV with different concentrations showed certain inhibition rates on BCPAP cells after 24, 48, and 72 hours of treatment, respectively. The IC 50 values of quercetin at 24, 48, and 72 hours were 0.361, 0.302, 0.291 µmol/ mL respectively (Fig. 6A). The IC 50 values of luteolin at 24, 48, and 72 hours were 0.425, 0.383, and 0.312 µmol/ mL respectively (Fig. 6B). The IC 50 values of β-sitosterol at 24, 48, and 72 hours were 0.361, 0.323, and 0.297 µmol/ mL respectively (Fig. 6C). The IC 50 values of kaempferol at 24, 48, and 72 hours were 0.311, 0.285, 0.262 µmol/mL respectively (Fig. 6D). Among them, kaempferol showed the strongest inhibition on the growth of BCPAP cells (IC 50 value was the lowest). Therefore, kaempferol was selected as the most effective compound of PV in subsequent experiments. It was revealed that kaempferol has a significant time effect relationship and a dose effect relationship on the inhibition of BCPAP cells, and the inhibition rate is the highest.

Discussion
PTC is the most common endocrine system tumor, which belongs to the category of "goiter and tumor" disease in TCM which accounts for about 80% to 90% of DTC. [14,15] TCM has shown great potential in the prevention and treatment of tumors and other diseases. The separation and extraction of effective components from TCM is an important means of new drug development. [16] It is reported that more than 80% of small molecule antitumor drugs are natural products and their derivatives. [17] The diverse pharmacological activities of natural compounds provide an important basis for the study of the mechanism of biological functions. Previous studies have shown that PV may inhibit the proliferation of tumor cells, induce apoptosis, and regulate immunity in many kinds of tumors. [12] However, due to the complexity of the components of TCM, its main effective components and related mechanisms are still unclear. Network pharmacology is to use systems biology, high-throughput screening and other technologies to systematically analyze the interaction network of "TCM components and disease targets," and to reveal the intervention and impact of TCM on body diseases as a whole, which is consistent with the theory of TCM treatment. [18] In the present study, the main active components, targets and functional pathways of PV in the treatment of PTC were preliminarily explored by network pharmacology. The results showed that there were 11 main active ingredients of PV in the treatment of PTC, including β-sitosterol, kaempferol, spinasterol, stigmasterol, delphinidin, luteolin, vulgaxanthin-I, poriferasterol monoglucoside_qt, stigmast-7-enol, morin, quercetin. According to the analysis of the network of the "active ingredient-target-disease" with Cyto scape 3.7.2, the top 4 active compounds of PV are quercetin, luteolin, β-sitosterol, kaempferol. Quercetin is a unique bioactive flavone and antioxidant, which plays a full role in reducing different human cancers that can directly promote apoptosis of tumor cells, so it can inhibit the progression of many human cancers. [19] It has been shown that quercetin can inhibit thyroid cancer cells in vitro. [20] In addition, it was also reported that quercetin is effective in the treatment of both medullary and papillary of human thyroid cancer. [21] β-sitosterol is the most abundant and widely distributed in vegetable foods rich in oil. It was found that prominent in vitro antiproliferative and pro-apoptotic effects of β-sitosterol in MDA-MB231 cells. [22] Previous study revealed that the silver nanoparticles of β-sitosterol can effectively induce the toxicity and early apoptosis of human colon cancer cells by enhancing the expression of p53 protein. [23] Luteolin may show effective cytotoxicity to thyroid cancer cells in vitro and in vivo by blocking BANCR/ TSHR signals which may become a potential important anticancer agent for thyroid cancer progression. [24] Meanwhile, luteolin exhibits an antitumor effect in the HCC cell line SMMC-7721 by suppress cell viability, induce G0/G1-phase arrest, and increase cellular apoptosis. [25] Kaempferol may inhibited cellular growth in a dose dependent manner in F9 and thyroid cancer cells. [26] It was reported that kaempferol inhibits the invasion and migration of renal cancer cells through the downregulation of AKT and FAK pathways. [27] Taken together, the 4 core components of PV showed different levels of anticancer activity.
These active ingredients act on 83 target proteins of PTC. PPI network shows that the genes with the top 6 connectivity are VEGFA, TP53, JUN, PTGS2, IL6, and IL1B, respectively. It shows that these targets are related to the proliferation and differentiation of tumor cells and tumor progression in the present study. Among them, VEGFA is a key mediator of cancer-associated neo-angiogenesis and progression which have been explored extensively for cancer therapy. It was found that the importance of increased VEGFA protein in tumors and may explain the higher protein quantification observed in the PTC tumor samples in comparison to the goiter. [28] It was also revealed that the molecular state of VEGFA may play an important role in the progression of PTC. VEGFA mRNA is overexpressed in thyroid carcinoma, especially in thyroid carcinoma with lymph node metastasis. [29] Previous research results found that the role of TP53-rs1042522 polymorphism in development of thyroid carcinoma which related to decreased risk of PTC, smaller tumor size, and lower incidence of vascular invasion. [30] In addition, it was indicated that TP53 may lead to the risk of PTC in individuals exposed to radiation in late childhood, adolescence or adulthood in previous study. [31] It has been demonstrated that the gene JUN was closely connected with PTC genesis which may benefit the cure of PTC patients. [32] It was found that the PTGS2 expression was significantly higher in patients with lymph node metastasis or extrathyroidal extension. [33] IL6 plays an important role in the progression of thyroid cancer that targeting IL6 signal transduction may help the clinical treatment of patients with thyroid cancer with more invasive tumor characteristics. [34] The IL6/JAK2/STAT3 pathway was considered to be associated with PTC malignant behaviors. [35] IL1B polymorphism may be associated with PTC risk and a predictor of lymph node metastasis of PTC patients in Korean population. [36] The results of molecular docking showed that the binding strength of key components corresponding to key targets was good. It further verified that PV may play a role in the treatment of PTC through these core pharmacodynamic components acting on the corresponding targets.
KEGG and GO analysis showed that the potential mechanism of PV in treating PTC mainly focused on cell proliferation, differentiation and anti-inflammatory immunity. The occurrence, development, malignant transformation, invasion, and metastasis of tumor are accompanied by inflammatory reaction. Inflammatory cells and molecules play a key role in the formation and maintenance of tumor immunity. [37] It had been demonstrated that the level of IL6 in patients with thyroid cancer is elevated, and there is a close relationship between thyroid cancer and inflammation. [38] The level of serum interleukin fully reflects this relationship and the potential mechanism of Table 3 The results of molecular docking.    factor kappa B, thus activating the expression of effector genes IL6 and vascular endothelial growth factor which promoting the occurrence and development of tumors. [39] PI3K-AKT signal pathway is involved in cell proliferation, apoptosis, cycle regulation and other pathophysiological processes. [40] P53 is one of the main apoptotic signal pathways, which can promote apoptosis by interacting with BCL-2 family proteins in the cytoplasm. [41] PV may treat PTC through multiple pathways and multiple targets. Besides above genes and pathways, it was also reported that BRAF gene mutation [42] and TERT promoter mutation [15] were closely associated with PTC recurrence and prognosis which should focus on the key points research in the future.

Conclusion
Eleven components and 83 corresponding targets were obtained in the component target network of PV, of which 6 were the core targets of PV in the treatment of PTC. Quercetin, luteolin, β-sitosterol, kaempferol may be the core components of PV in the treatment of PTC. VEGFA, TP53, JUN, PTGS2, IL6, and IL-1B were considered to be the targets for the treatment of PTC. IL-17 signaling pathway and PI3K-Akt signaling pathway may affect the recurrence and metastasis of PTC. In conclusion, this study used network pharmacology and molecular docking technology to clarify that PV may play a role in the treatment of PTC through multiple components, multiple targets and multiple pathways. In addition, further vivo experiment was carried to verify it. Kaempferol is consider to be the most useful compound which may reduce the protein expression levels of IL6, VEGFA, JUN, TP53, 1L-1B, and PTGS2, respectively. This study laid a foundation for the in-depth study of the material basis and mechanism of action of PV, and provided a scientific basis for the clinical application of PV.